**Traffic and Marijuana Policy Bundle
**July 5th, 2023
*scale correction

clear all

import excel "C:\Users\parks\OneDrive - The Pennsylvania State University\Desktop\Cannabis Legalization and Policy Bundles Projects\Traffic Death\Cannabis bundles LR _mgpark.xlsx", sheet("Data") firstrow

**variable
*depedent
*Fatality rate per 10 Billion Vehicle Miles Traveled
gen ftl_rate=Fatality/TotalMiles*10000

*Drunken driver crash rate per 10 Billion Vehicle Miles Traveled
gen alcoholcrash_rate=alcoholcrashes/TotalMiles*10000

*Pedestrian fatality rate per 10 Billion Vehicle Miles Traveled
gen pedestrian_rate=Pedestrian/TotalMiles*10000

*Total occupant fatality rate per 10 Billion Vehicle Miles Traveled
gen occ_rate=TotalOccupantsKilled/TotalMiles*10000

*Passenger cars fatality rate per 10 Billion Vehicle Miles Traveled
gen occ_car_rate=PassengerCars/TotalMiles*10000

*Light trucks fatality rate per 10 Billion Vehicle Miles Traveled
gen occ_litruck_rate=LightTrucks/TotalMiles*10000

*Motorcycles fatality rate per 10 Billion Vehicle Miles Traveled
gen occ_motor_rate=Motorcycles/TotalMiles*10000


*population controlled
*Fatality rate per 10 Million population
gen ftl_rate2=Fatality/Population*10000000

*Drunken driver crash rate per 10 Million population
gen alcoholcrash_rate2=alcoholcrashes/Population*10000000

*Pedestrian fatality rate per 10 Million population
gen pedestrian_rate2=Pedestrian/Population*10000000

*Total occupant fatality rate per 10 Million population
gen occ_rate2=TotalOccupantsKilled/Population*10000000

*Passenger cars fatality rate per 10 Million population
gen occ_car_rate2=PassengerCars/Population*10000000

*Light trucks fatality rate pper 10 Million population
gen occ_litruck_rate2=LightTrucks/Population*10000000

*Motorcycles fatality rate per 10 Million population
gen occ_motor_rate2=Motorcycles/Population*10000000


*indepedent
*Alcohol consumption – NIAAA
gen alcoh=AlcoholPercapitaconsumptionA

*Income per capita in thousands of dollars
gen income=Personalincomepercapita/1000

*Unemployment rate
gen unemploy=Unemploymentrate

*Population density
gen popdensity=Populationdensity

*Temperature – NOAA
gen temp=Temperature

*Precipitation - NOAA
gen preci=Precipitation


*policy variables
*08BAC .08 BAC laws limit the number of drinks a driver 
gen bac=AU 
gen bacfrac=AV 

*Open container law 
gen opencontainer=OpenCont 

*Seat belt laws
gen seat_prim=beltprim 
gen seat_second=beltscnd

**descriptive
su ftl_rate ftl_rate2 alcoholcrash_rate alcoholcrash_rate2 pedestrian_rate pedestrian_rate2 occ_rate occ_rate2 occ_car_rate occ_car_rate2 occ_litruck_rate occ_litruck_rate2 occ_motor_rate occ_motor_rate2 PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim
ta yearcode

**heteroskedasticity detection
qui: reg ftl_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year i.FIPS
predict residuals1, residuals
hettest residuals1

qui: reg alcoholcrash_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year i.FIPS
predict residuals2, residuals
hettest residuals2

qui: reg pedestrian_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year i.FIPS
predict residuals3, residuals
hettest residuals3

qui: reg occ_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year i.FIPS
predict residuals4, residuals
hettest residuals4

qui: reg occ_car_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year i.FIPS
predict residuals5, residuals
hettest residuals5

qui: reg occ_litruck_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year i.FIPS
predict residuals6, residuals
hettest residuals6

**regression
xtset FIPS yearcode

xtgls ftl_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)
xtgls alcoholcrash_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)
xtgls pedestrian_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)
xtgls occ_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)
xtgls occ_car_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)
xtgls occ_litruck_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)

qui: xtgls ftl_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)
est sto m1
qui: xtgls alcoholcrash_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)
est sto m2
qui: xtgls pedestrian_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)
est sto m3
qui: xtgls occ_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)
est sto m4
qui: xtgls occ_car_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)
est sto m5
qui: xtgls occ_litruck_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)
est sto m6
esttab m1 m2 m3 m4 m5 m6, cells(b(fmt(3)) ci(fmt(3) par) p(fmt(3) par)) scalars(ll_0 ll chi2) aic bic ///
gaps lines collabels(none) ///
addnotes("95% confidence intervals in brackets; p-values in parentheses.")

**multiple correction wrong false
*The Bonferroni correction is a conservative method to control for multiple testing, reducing the risk of false positives
qui: xtgls ftl_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)
test PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim, mtest(bonferroni)
lincom PharmImpAnyScale
lincom PermImpAnyScale
lincom FiscalImpAnyScale
lincom alcoh
lincom income 
lincom unemploy 
lincom popdensity 
lincom temp 
lincom preci 
lincom bac 
lincom opencontainer 
lincom seat_prim

qui: xtgls alcoholcrash_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)
test PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim, mtest(bonferroni)

qui: xtgls pedestrian_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)
test PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim, mtest(bonferroni)

qui: xtgls occ_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)
test PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim, mtest(bonferroni)

qui: xtgls occ_car_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)
test PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim, mtest(bonferroni)

qui: xtgls occ_litruck_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)
test PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim, mtest(bonferroni)

**multi collinearity
qui: reg ftl_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim
estat vif

qui: reg alcoholcrash_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim
estat vif

qui: reg pedestrian_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim
estat vif

qui: reg occ_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim
estat vif

qui: reg occ_car_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim
estat vif

qui: reg occ_litruck_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim
estat vif

*aic bic
qui: xtgls ftl_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year
est sto m1
qui: xtgls alcoholcrash_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year
est sto m2
qui: xtgls pedestrian_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year
est sto m3
qui: xtgls occ_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year
est sto m4
qui: xtgls occ_car_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year
est sto m5
qui: xtgls occ_litruck_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year
est sto m6
esttab m1 m2 m3 m4 m5 m6, cells(b(fmt(3)) ci(fmt(3) par) p(fmt(3) par)) scalars(ll_0 ll chi2) aic bic r2 ///
gaps lines collabels(none) ///
addnotes("95% confidence intervals in brackets; p-values in parentheses.")

**post estimation
qui: xtgls ftl_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)
margins, at(PharmImpAnyScale=(0(8.3)99.6))
marginsplot, xtitle("Pharmaceutical scale") graphregion(fcolor(white))
graph export 1.png, replace

margins, at(PermImpAnyScale=(0(8.3)99.6))
marginsplot, xtitle("Permissive scale") graphregion(fcolor(white))
graph export 2.png, replace

margins, at(FiscalImpAnyScale=(0(8.3)99.6))
marginsplot, xtitle("Fiscal scale") graphregion(fcolor(white))
graph export 3.png, replace


qui: xtgls alcoholcrash_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)
margins, at(PharmImpAnyScale=(0(8.3)99.6))
marginsplot, xtitle("Pharmaceutical scale") graphregion(fcolor(white))
graph export 4.png, replace

margins, at(PermImpAnyScale=(0(8.3)99.6))
marginsplot, xtitle("Permissive scale") graphregion(fcolor(white))
graph export 5.png, replace

margins, at(FiscalImpAnyScale=(0(8.3)99.6))
marginsplot, xtitle("Fiscal scale") graphregion(fcolor(white))
graph export 6.png, replace


qui: xtgls pedestrian_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)
margins, at(PharmImpAnyScale=(0(8.3)99.6))
marginsplot, xtitle("Pharmaceutical scale") graphregion(fcolor(white))
graph export 7.png, replace

margins, at(PermImpAnyScale=(0(8.3)99.6))
marginsplot, xtitle("Permissive scale") graphregion(fcolor(white))
graph export 8.png, replace

margins, at(FiscalImpAnyScale=(0(8.3)99.6))
marginsplot, xtitle("Fiscal scale") graphregion(fcolor(white))
graph export 9.png, replace


qui: xtgls occ_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)
margins, at(PharmImpAnyScale=(0(8.3)99.6))
marginsplot, xtitle("Pharmaceutical scale") graphregion(fcolor(white))
graph export 10.png, replace

margins, at(PermImpAnyScale=(0(8.3)99.6))
marginsplot, xtitle("Permissive scale") graphregion(fcolor(white))
graph export 11.png, replace

margins, at(FiscalImpAnyScale=(0(8.3)99.6))
marginsplot, xtitle("Fiscal scale") graphregion(fcolor(white))
graph export 12.png, replace


qui: xtgls occ_car_rate PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)
margins, at(PharmImpAnyScale=(0(8.3)99.6))
marginsplot, xtitle("Pharmaceutical scale") graphregion(fcolor(white))
graph export 13.png, replace

margins, at(PermImpAnyScale=(0(8.3)99.6))
marginsplot, xtitle("Permissive scale") graphregion(fcolor(white))
graph export 14.png, replace

margins, at(FiscalImpAnyScale=(0(8.3)99.6))
marginsplot, xtitle("Fiscal scale") graphregion(fcolor(white))
graph export 15.png, replace

qui: xtgls occ_litruck_rate  PharmImpAnyScale PermImpAnyScale FiscalImpAnyScale alcoh income unemploy popdensity temp preci bac opencontainer seat_prim i.year, panels(correlated)
margins, at(PharmImpAnyScale=(0(8.3)99.6))
marginsplot, xtitle("Pharmaceutical scale") graphregion(fcolor(white))
graph export 16.png, replace

margins, at(PermImpAnyScale=(0(8.3)99.6))
marginsplot, xtitle("Permissive scale") graphregion(fcolor(white))
graph export 17.png, replace

margins, at(FiscalImpAnyScale=(0(8.3)99.6))
marginsplot, xtitle("Fiscal scale") graphregion(fcolor(white))
graph export 18.png, replace
